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Multiple-time scaling and universal behaviour of the earthquake inter-event time 

distribution 
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The inter-event time distribution characterizes the temporal occurrence in seismic catalogs. Uni- 
versal scaling properties of this distribution have been evidenced for entire catalogs and seismic 
sequences. Recently, these universal features have been questioned and some criticisms have been 
raised. We investigate the existence of universal scaling properties by analysing a Californian catalog 
and by means of numerical simulations of an epidemic-type model. We show that the inter-event 
time distribution exhibits a universal behaviour over the entire temporal range if four characteristic 
times are taken into account. The above analysis allows to identify the scaling form leading to 
universal behaviour and explains the observed deviations. Furthermore it provides a tool to identify 
the dependence on mainshock magnitude of the c parameter that fixes the onset of the power law 
decay in the Omori law. 
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Seismic occurrence is a phenomenon of great complex- 
ity involving different processes acting on different time 
and space scales. In the last decade a unifying picture 
of seismic occurrence has been proposed via the inves- 
tigation of D(At), the distribution of inter-event times 
At between successive earthquakes These studies 

have shown that, rescaling inter-event times by the aver- 
age occurrence rate R, D (At) follows the scaling relation 



D(At) = Rf(RAt) 



(1) 



where the functional form of f(x) is quite independent of 
the geographic zone and the magnitude threshold. The 
above relation suggests that R is a non universal quan- 
tity and is the only typical inverse time scale affecting 
D(At). This result, obtained for periods of stationary 
rate, has been generalized to non stationary periods [4[ 
and Omori sequences [jj @ . The scaling relation JTJ has 
been also observed for volcanic earthquakes Q- On the 
other hand, recent studies have questioned the universal- 
ity of the inter-event time distribution @~ 14 1 . In partic- 



ular, deviations from universality at small At have been 
related to the interplay between correlated earthquakes, 
following a gamma distribution, and uncorrelated events, 
following a pure exponential decay [l4| . This behaviour is 
well described by numerical simulations of the Epidemic 
Type Aftershock Sequence (ETAS) model [15j. Indeed, 
analytical studies [9| and a previous numerical analysis 
of the ETAS model Q, have shown that the functional 
form of D(At) depends on the ratio between correlated 
and independent earthquakes /C. The problem has been 
also attacked within the theoretical framework of prob- 
ability generating functions 1 lM 1 31] . Saichev & Sornette 



in Eq.|T]) depends on K, and on some other parameters 
of the model. This behavior is confirmed if non-linear 
contributions are taken into account (l3j ] . 

Multiplicity of characteristic times is often observed in 
the dynamics of complex systems, where different tem- 
poral scales are associated to the relaxation of different 
spatial regions or structures. For instance, their exis- 
tence is a well established property in glassy materials, 
polymers or gelling systems, where they originate from 
the relaxation of complex structures at different meso- 
scopic scales, or else from the emergence of competing 
interactions 16[. Moreover, the coexistence of different 



(SS) have obtained an exact nonlinear integral equation 
for the ETAS model and solved it analytically at linear 
order [ll, 12 1. This solution shows that the function f(x) 



physical mechanisms acting at different spatio-temporal 
scales may also give rise to complex temporal scaling [ijj ■ 
Therefore, the identification of the number of relevant 
time scales controlling universal behaviours is a very de- 
bated subject in complex systems. 

In this paper, we do not assume the existence of a 
unique time scale 1/R, as in Eq.([T]). We show that 
four typical timescales are relevant for the inter-event 
time distribution scaling: the inverse rate of independent 
events /i, the average inverse rate of correlated events, the 
time parameter c defined in the Omori law and the cata- 
log duration T. These different time scales lead to devia- 
tions from the simple scaling ([T]). Nevertheless, we show 
that the inter-event time distribution can be expressed in 
a universal scaling form in terms of these four character- 
istic times. This scaling form allows to better enlighten 
the mechanism leading to universality for D(At) and the 
deviations from it. The above analysis also clarifies the 
dependence of c on the mainshock magnitude for inter- 
mediate mainshock sizes. 

We assume that seismic occurrence can be modeled 
by a time-dependent Poisson process with instantaneous 
rate A(t). In this case, the inter-event time distribution 
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for the temporal interval [0, T] is [B[ 



D(At) 
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■ T-At 



dsX(s)X(s + A<)exp - ^ 



>+At 



+ A(At)exp-/o At£inA ^ 



duX(u) 



(2) 



where N = J Q dsX(s) is the number of events. The 
widely accepted scenario is that seismic occurrence can 
be considered as the superposition of Poissonian events 
occurring at constant rate fi and independent aftershock 
sequences, which gives 



A(i) = /i+ J2 ViiP- 1 ) 



i:ti<t 



(3) 



where p > 1 is the exponent of the Omori law. The quan- 
tity cii is proportional to the rate of aftershocks correlated 
to the i-th mainshock since, from Eq.©, the total num- 
ber of events triggered by the i-th mainshock is cijCj. The 
productivity law [l8| indicates that at is exponentially 
related to the mainshock magnitude <Zj = AlO anii . The 
ETAS model assumes Cj = c whereas recent studies on 
experimental catalogs have obtained a — clO ami which 
leads to the so-called generalized Omori law ■ 
The dependence of on mainshock magnitudes has been 
attributed to a dynamical scaling relation involving time, 
space and energy [2oM23l | or to catalog incompleteness 
24}. Inserting Eq.© in Eq.©, we obtain a scaling form 
for D(At) expressing time in unit of l//x, 



D(At) = fxG (/j,At, fi/a, fic, fiT) 



(4) 



where a (c) is the value of a,i (ci) averaged over all main- 
shocks. We show that Eq.© represents a particular case 
of the more general scaling form ©. Indeed, by defini- 
tion, R in the time interval [0, T] is the inverse of the av- 
erage At, and from Eq.© R — fi/H (fi/a, [ic, /xT), with 
H(y,z,w) — J dxG(x,y, z,w)x. Therefore, expressing /i 
in terms of R and setting K. — a c, we obtain 



D(At) = Rd (RAt, K,, nc, /jT) . 



(5) 



For the ETAS model JC = He is the branching ratio, i.e. 
the number of direct aftershocks per earthquake. 

The complex form of Eq. ((3]) does not allow the 
full derivation of an analytical expression for D(At), 
unless one uses specific assumptions. SS [12j, for 
instance, have exactly calculated D(At) for T —> 
oo in the hypothesis that each earthquake triggers, 
on average, the same number of aftershocks, i.e. 
a,i = A and Ci = c. This solution exactly fol- 
lows the scaling form Eq.© with Gi(x, y, z, oo) = 
exp -y)- {uj 2 -p - 1) zy/((l - y)(2 -p))) 

((l - y + yuj 1 ^) 2 + (p - l)uj- p y{l - y) / 'zj withw = 



1 + (1 — y)x/z. For a given choice of the parameters 
/C, \ic and p the above expression leads to a D(At) in 
good agreement with the experimental distribution. In- 
terestingly, the above expression coincides with the linear 
order of the ETAS model expansion, in the limit At 3> c. 
Higher order terms lead to small differences with the 
above solution (l3j . 

A useful example to understand the role of the dif- 
ferent time scales in D(At) can be obtained if we limit 
the calculation to events in a single aftershock sequence. 
A scaling form consistent with Eq.© has been already 
obtained in ref.[H assuming (i = 0. Here, we restrict to 
large T and At <C T assuming that A(t) is about constant 
in At. This choice does not represent a loss of general- 
ity for sufficiently small At. Under these assumptions, 
Eq.© becomes D(At) = i J Q T ds\ 2 (s) exp~ x ^ At which 
can be analytically integrated. Three typical regimes can 
be identified: i) At large times fiAt > 1 (i.e. x > l),with 
At/T <C 1, D(At) decays as exp (— fiAt) as already ob- 
served [9j. ii) At intermediate times c < At < 1//J, 
(y < x < 1), we observe a power law decay At~ e with 
e = 2 — 1/p, as predicted by Utsu 25] . iii) At small times 
At <C c (x/y <C 1), D(At) becomes At independent. 
The three regimes can be identified in Fig. [T]where, for a 
fixed value of fx and different c, we plot AtD(At) vs RAt. 
This is equivalent to the representation adopted by [14 1 
and allows to better enlighten deviations from the scaling 
relation ©. We observe that all curves present a peak 
at RAt ~ 1 and then exponentially decay for RAt > 1. 
Conversely, at small times (At < c), all the curves in- 
crease linearly since D(At) is constant. The intermediate 
regime can be observed only for the two smallest values 
of c, since only in these cases /ic -C 1 and the intermedi- 
ate regime has a finite extension. In this regime an about 
flat behavior is observed since for p = 1.05, e ~ 1. Notice 
that one of the theoretical curves provides a good quali- 
tative fit of the distribution obtained from experimental 
data 
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In Fig. 2 we present the results of numerical simula- 
tions of the ETAS model obtained following the method 
of ref . |27| . We perform extensive simulations in order to 
recover the limit T — > oo and neglect the dependence 
on (j,T. Previous studies [13] have proposed that D(At) 
only depends on the branching ratio and used this result 
to measure the ratio between triggered and independent 
events in experimental catalogs. Sornette et al. [13| have 
shown that the functional form of D(At) also depends on 
p, whose experimental values fluctuate around one. Our 
numerical simulations with a random p g]1, 1.6] indicate 
that D(At) depends on the average p value. In the fol- 
lowing we present result for simulations with p = 1.2. 
A very similar pattern is obtained for other p values. 
Fig. [2^ shows that for fixed values of fi and c the curves 
exhibit different behaviours for different K., in agreement 
with previous results. We then focus on the role of the 
parameter c. In Fig. [2b, we plot AtD(At) at constant 
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FIG. 1: AtD(At) vs. RAt for the single Omori sequence 
with fi = 3 • 10 -7 events per second and varying c from 40 to 
4 ■ 10 7 sec (from left to right). Other model parameters are 
T = 1 • 10 7 sec, A — 1 and p = 1.05. Open circles represent 
AtD(At) for the California catalog [261 ]. 



JC and c and for different values of //, as in ref.|l4|. We 
confirm the existence of deviations from scaling ([1} at 
small At which can be attributed to c. In order to show 
that the dependence on c enters in the scaling form as 
/ic, in Fig. [2t, we fix the branching ratio IC and vary /j, 
as in Fig. [2b, but now c is allowed to vary keeping con- 
stant lie. In this case all the distributions collapse on the 
same master curve revealing an "universal" behaviour 
also at small At. This result indicates that deviations 
from Eq.(fT]) rely on the presence of the variable /ic in 
Eq. (J4j . The dependence on this variable is relevant for 
At < c and becomes negligible at larger At. This ac- 
counts for the appearance of deviations from universality 
only at small At. 

We now explore the role of the characteristic time c 
on the scaling properties of AtD(At) for the experimen- 
tal catalog [26(. We first consider the whole catalog. 
In this case [iT is very large and does not affect the 
scaling form ([5]). To isolate the dependence on /ic, we 
plot AtD(At) including in the analysis only earthquakes 
above a lower magnitude thresholds m t h- p/a, indeed, 
should not depend on m t h [HI and AtD(At) is expected 
to depend only on RAt and flc. Fig. 3a shows deviations 
from universality. Since these deviations are confined at 
small At they are not easy to detect in the usual plot 
R-^DiAt) vs RAt 0. Deviations must be attributed 
to fie and allow to identify c from the crossover points 
separating the linear growth from the plateau (identi- 
fied by arrows in Fig. 3a). Using the known values of 
R we obtain that c is quite independent of m t h- This 
is consistent with a independent of mainshock magni- 
tudes but also with c t oc I0 a ( m *-™ th ) and a < b pjlj . 
The dependence of Cj on and its influence on the 
D(At) can be obtained by restricting the analysis to 
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FIG. 2: AtD(At) vs. RAt for the simulated catalogs with 
T = 10 7 and p = 1.2. (a) AtD(At) for c = 3000, ju = 8 • 10~ 6 
events per unit time and varying K, using different A values 
reported in the legend, (b) AtD(At) for c = 3000, K = 0.5 
and varying n according to the legend in panel (c). (c) We fix 
c/j, = 0.24 obtaining data collapse with K. = 0.5. In all figures 
the time unit is the iteration step of the model. 



temporal periods soon after mainshocks. In these inter- 
vals, /j, can be neglected simplifying the scaling relation 
(0 . We further observe that for a single Omori sequence 
R = N/T = a lCl [l - (T/c,;) (1 - p) ] /[T(p - 1)] and there- 
fore R/a,i is a function only of the ratio T/ci. In this case 
the scaling simplifies to 



D(At) = RG 2 \ RAt, ^ 



(6) 



We start by considering the main-aftershock sequences 
for the three largest shocks recorded in the catalog: Lan- 
ders, Northridge and Hector Mine. We consider as after- 
shocks all events with m > 2.5 occurring in a temporal 
window [0,T] after the main event and within the after- 
shock zone, i.e. a radius L = 0.01 • 10 O 5m4 km from the 
mainshock. Different definitions of the aftershock zone 
[29j | lead to very similar results. We first fix T — 10 days 
for all sequences. In Fig. [3b the curves do not collapse 
but show a progressive shift as the mainshock magnitude 
increases. This effect can be attributed to dependence of 
Cj on m,. Therefore, at fixed T the variable T/ci assumes 
different values for each sequence violating the collapse 
Eq.(fT]). As an alternative approach we use the criterion 
proposed in ref. [f| to identify the end of a sequence: 
namely a sequence ends when the rate X(t) reaches the 
average Poisson rate fi ~ 2 events/day. Fig. [3t clearly 
indicates a very good data collapse for the different se- 
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FIG. 3: (a) AID (At) vs. RAt for different magnitude thresh- 
olds mth; (b) for the Landers, Northridge and Hector Mine 
sequences (mainshock magnitudes m = 7.3,6.7,7.1, respec- 
tively) with fixed duration T = 10 days; (c) for the Landers, 
Northridge and Hector Mine sequences with T chosen follow- 
ing the criterion of ref.Q; (d) for all sequences with a main- 
shock magnitude m greater than 4 grouped in classes of m, 
m e [M, M + 0.5[. 



quences in good agreement with Eq.([T]). This can be un- 
derstood from the scaling relation © where the variable 
T/ci assumes the same value for all sequences. Indeed, 
according to the Omori law, for (>c the occurrence rate 
can be expressed as A(t) ~ lO Qmi t~z> and the condition 
A(T) = 2 provides T ~ 10 Qm */f>. For the largest main- 
shocks, Cj ~ I0 bnii /p [5], [|| and therefore, the condition 
a ~ b/p implies that the ratio T/ci is almost constant 
for the three sequences. 

Next, we extend the above analysis to all sequences 
with mainshock magnitude m > 4. To improve the 
statistics, we group mainshocks in classes of magnitude 
in 6 [M,M + 0.5[. Mainshocks are identified with the 
criterion suggested in ref.Q and the duration i s ag ain 
fixed by the condition //(T) = 2. Other methods 29] for 
aftershock identification provide similar results. Fig. [3pl 
shows data collapse for all M values. Since the crite- 
rion n(T) — 2 roughly implies T ~ 10 qM , the collapse 
of Fig. [5U, suggests that the dependence of Cj on the 
mainshock magnitude as a ~ I0 bmi /p is valid also for 
intermediate mainshock magnitudes. 

In conclusion, we address recent criticisms to the uni- 
versal behaviour of the intcr-evcnt distribution. We fol- 



low the approach of ref.[14j and show that D(At) does 
exhibit universal features on the whole temporal range if 
four characteristic time scales are taken into account. In 
particular, deviations at small At can be attributed to c 
scaling differently from 11. Whereas, by keeping constant 
T/c for different sequences, the D(At) collapse onto a 
unique master curve. 
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